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Heavy quarkonia are observed to be suppressed in relativistic heavy ion collisions relative to their 
production in p + p collisions scaled by the number of binary collisions. In order to determine if this 
suppression is related to color screening of these states in the produced medium, one needs to account 
for other nuclear modifications including those in cold nuclear matter. In this paper, we present new 
measurements from the PHENIX 2007 data set of J/tp yields at forward rapidity (1.2 < \y\ < 2.2) 
in Au+Au collisions at ^/sj^ = 200 GeV. The data confirm the earlier finding that the suppression 
of J/ip at forward rapidity is stronger than at midrapidity, while also extending the measurement to 
finer bins in collision centrality and higher transverse momentum (pt). We compare the experimental 
data to the most recent theoretical calculations that incorporate a variety of physics mechanisms 
including gluon saturation, gluon shadowing, initial-state parton energy loss, cold nuclear matter 
breakup, color screening, and charm recombination. We find J/V' suppression beyond cold- nuclear- 
matter effects. However, the current level of disagreement between models and d-|-Au data precludes 
using these models to quantify the hot-nuclear-matter suppression. 



PACS numbers: 25.75.Dw 

Heavy quarkonia have long been proposed as a sensi- 
tive probe of the color screening length and deconfine- 
ment in the quark-gluon plasma The picture that 
was originally proposed is complicated by other compet- 
ing effects that modify quarkonia production and sur- 
vival in cold and hot nuclear matter. The large sup- 
pression of J/tp in Pb-I-Pb collisions at ^/s^ — 17.2 
GeV measured by the NA50 experiment indicated sup- 
pression beyond that projected from cold nuclear matter 
effects and led to the initial conclusion that color screen- 
ing was the dominant mechanism • The expectation 
was that at the Relativistic Heavy Ion Collider (RHIC), 
where higher energy densities and temperatures are cre- 
ated, the J/tP suppression would be stronger and turn 
on in even more peripheral collisions. However, mea- 
surements from the PHENIX experiment's 2004 data set 
in Au-I-Au collisions at ^s„„ = 200 GeV revealed that 
the centrality-dependent nuclear modification factor Raa 
at midrapidity was the same within statistical and sys- 
tematic uncertainties as the NA50 result [B-Q- In addi- 
tion, for more central collisions (A^part ^ 100) the sup- 
pression was measured to be larger at forward rapidity 
(1.2 < |y| < 2.2) compared with midrapidity (Sj. This is 
opposite to the expectation that the suppression should 
be less at forward rapidity, where the energy density is 
lower. 

The initial estimates of cold nuclear matter (CNM) 
effects in the NA50 Pb+Pb data were based on p-\-A 
measurements [8] at higher collision energies. More re- 
cently, measurements of J/iJj production in p-f A colli- 
sions at ^/s^ =17.2 GeV made by NA60 have shown 
that CNM effects are stronger than the initial estimates. 
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due to a substantially larger J/ip effective breakup cross 
section Q . This resulted in a reduction of the estimated 
suppression due to possible hot nuclear matter effects 
from ^50% to ^-^25%, relative to cold nuclear matter ef- 
fects. Based on a systematic study of the energy depen- 
dence [l3] the breakup cross section is expected to be 
much smaller at the higher RHIC energies, leading to 
smaller overall CNM effects. Thus it now seems that at 
RHIC the suppression beyond cold nuclear matter effects 
at midrapidity could be higher compared to the NA50 
results, as one would expect due to the higher energy 
density. However, the question of why the observed sup- 
pression at RHIC is stronger at forward rapidity than at 
midrapidity remains less understood. 

Since that first PHENIX measurement, many alterna- 
tive explanations have been proposed which now require 
rigorous confrontation with the full set of experimental 
measurements. In this paper, we detail the measure- 
ment of forward rapidity J/^ yields and modifications in 
Au-I-Au colHsions at y^s^ — 200 GeV from data taken 
by the PHENIX experiment in 2007. The data set is 
more than three times larger than the earlier published 
2004 data set. In addition, significant improvements in 
our understanding of the detector performance and signal 
extraction have led to a reduction in systematic uncer- 
tainties. We present details of this new analysis as well 
as comparisons with theoretical calculations that have 
attempted to reconcile the earlier data in terms of com- 
peting physics mechanisms. 



I. DATA ANALYSIS 

The PHENIX experiment is described in detail in [ij] . 
For the forward rapidity J/tp data analysis presented 
here, the PHENIX experiment utilizes one global detec- 
tor, the beam-beam counter (BBC), for event central- 
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ity characterization and z-vertex determination, and two 
muon spectrometers (North and South) for measuring 
the J 14' via the dimuon decay channel The BBC is de- 
scribed in detail in fT?]. It comprises two arrays of 64 
quartz Cerenkov counters that measure charged parti- 
cles within the pseudorapidity range (3.0 < \rj\ < 3.9). 
The BBC is also used as the primary Level-1 trigger 
for Au-|-Au minimum bias events. The two muon spec- 
trometers are described in detail in , and comprise an 
initial hadronic absorber followed by three sets of cath- 
ode strip chambers in a magnetic field, referred to as the 
Muon Tracker (MuTR). Finally, there are five planes of 
active larocci tubes interleaved with additional steel ab- 
sorber plates, referred to as the Muon Identifier (MuID). 
The muon spectrometers measure J/ tps over the rapidity 
range 1.2 < |?;| < 2.2. 

The PHENIX data acquisition system is capable of 
recording minimum bias Au-|-Au collisions at high rates 
(> 5 kHz) with a data archival rate in excess of 600 
MB/s. During the 2007 Au+Au run at ^/i;;^ = 200 
GeV, the experiment recorded 82% of all collisions where 
the minimum bias level-1 trigger fired. Therefore no ad- 
ditional muon-specific trigger was necessary. After data 
quality cuts to remove runs where there were significant 
detector performance variations, we analyzed 3.6 x 10^ 
minimum bias Au-|-Au events. The BBC Level-2 trigger 
used as the minimum bias trigger for Au-|-Au events re- 
quired at least two hits in each of the BBC arrays and 
a fast-reconstructed z-vertex within ±30 cm of the nom- 
inal center of the detector. An additional selection for 
our minimum bias definition in offline reconstruction in- 
cludes a requirement of at least one neutron hit in each 
of our Zero Degree Calorimeters (ZDCs). This removes 
a 1-2% background event contamination in the periph- 
eral event sample. The minimum bias definition corre- 
sponds to 92% ± 3% of the inelastic Au-|-Au cross sec- 
tion 0. We further categorize the events in terms of 
centrality classes by comparing the combined North and 
South BBC charge to a negative binomial distribution 
for the number of produced particles within the BBC 
acceptance combined with a Glauber model of the colli- 
sion [iJl- For each centrality category, the mean num- 
ber of participating nucleons (A'part), binary collisions 
(TVcou), and impact parameter (&), and their associated 
systematic uncertainties are shown in Table ID 

From this minimum bias data sample, we recon- 
struct muon candidates by finding tracks that penetrate 
through all layers of the MuID, then matching to tracks 
in the MuTR. The requirement of the track penetrating 
the full absorber material through the MuID significantly 
reduces the hadron contribution. However, there is a 
small probability of order ~1/1000 for a charged hadron 
to penetrate the material without suffering a hadronic in- 
teraction. This is referred to as a punch-through hadron. 
Additionally, the current muon spectrometer cannot re- 
ject most muons that originate from charged pions and 
kaons which decay before the absorber in front of the 
MuTR. Pairs of muon candidate tracks are selected and 



TABLE I: Mean A'part, A^'coii, and impact parameter values 
and systematic uncertainties in each centrality category. 
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a combined fit is performed with the collision z-vertex 
from the BBC. We apply various cuts to enhance the 
sample of good muon track pairs, including cuts on the 
individual track values, the matching between position 
and direction vectors of the MuID track and the MuTR 
track projected to the front of the MuID, and finally the 
track pair and BBC z-vertex combined fit Xvtx- 

We then calculate the invariant mass of all muon can- 
didate unlike charged sign pairs in various bins in rapid- 
ity, pt, and centrality. Due to the high particle mul- 
tiplicity in Au-|-Au events, there is a significant back- 
ground under the J/ip peak. In the 0-5% most central 
Au-|-Au collisions, within the mass window around the 
J/V' (2.6 < M [GeV/c^] < 3.6) the signal to background 
is of order 2.9% in the South muon spectrometer and 
0.7% in the North muon spectrometer (which has a dif- 
ferent geometric acceptance and significantly higher oc- 
cupancy). The background comprises two components. 
First, there is the combinatorial background from uncor- 
relatcd track pairs. Second, there is a correlated back- 
ground from physical sources including open charm pair 
decays (e.g. + D" {K^ ^i+v^,) + (K+fi'T^), open 
beauty pair decays, and Drell-Yan. The combinatorial 
background is estimated and subtracted by using event 
mixing to recreate the background from uncorrelated 
pairs. Pairs are created from different Au-|-Au events 
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within the same category in Au+Au centrahty and BBC 
coUision z- vertex. The mixed-event invariant mass distri- 
butions are calculated for unlike-sign and like-sign pairs. 
We treat the real event like-sign pairs as being purely 
from combinatorial background, since the contribution 
of the above-mentioned correlated physical background 
is negligible relative to the combinatorial background. 
Thus, we determine the mixed event normalization for 
the unlike-sign case by calculating the normalization fac- 
tor between the mixed event and real event like-sign 
counts. We have confirmed that the mixed event and 
real event like-sign invariant mass distributions match 
over the full mass range used in the analysis. 

The total J/V' counts recorded in all Au-|-Au colli- 
sions are ^9100 and ~4900 in the south and north muon 
spectrometers, respectively. Two example unlike-sign in- 
variant mass distributions before (upper panel) and af- 
ter (lower panel) mixed-event combinatorial background 
subtraction are shown in Fig.[T] The number of measured 
J 1^ is derived from the subtracted spectra by fitting to 
the data the line shape of the J/^, as determined in 
p+p collisions [l5j, and an exponential for the remaining 
correlated physical background. Note that in the higher 
Pt > 2 GeV/c bins, there is some small acceptance for 
the p, Lo and </> and thus these additional components are 
included in the fit. At low pT, the acceptance goes to 
zero for lower invariant mass due to the minimum re- 
quired momentum for each muon to penetrate the MuID 
and the angular acceptance of the spectrometer. This 
is accounted for by folding the fit function with an ac- 
ceptance function that is calculated from Monte Carlo 
simulation and goes to zero at small invariant mass, as 
expected. We perform a set of fits where we vary the 
invariant mass range, the line shape of the J/V'j and the 
normalization of the mixed-event sample by ±2% to de- 
termine the mean extracted J/V' signal and systematic 
uncertainty from the RMS of the different results. Note 
that any bin where the extracted signal is of less than one 
standard deviation significance (including statistical and 
systematic uncertainties) is quoted as a 90% confidence 
level upper limit (CLUL) based on Poisson statistics. 

We also estimated the combinatorial background using 
the like-sign method, i.e. TV™"'' = 2^iV++ x iV__. In 
this case there is no event mixing, the background is esti- 
mated purely from same-event like-sign pairs (instead of 
mixed-event pairs). The two methods agree over most of 
the centrality range; however, for the more central events 
the like-sign method results in somewhat lower extracted 
counts {'^ 10%). We take the average of the two signal 
extraction methods and assign an additional systematic 
uncertainty due to the difference, although for peripheral 
bins it is a negligible difference. 

We then calculate the J/ip invariant yield for each cen- 
trality bin and also in bins in pT by the following equa- 
tion: 



where B^^ is the branching fraction of J /ip to muons, 
Nj/^ is the number of measured J/V', Nevt is the num- 
ber of events in the relevant Au+Au centrality cate- 
gory, Ae is the detector geometric acceptance and effi- 
ciency, and Apr and Ay are the bin width in pT and 
y. For the pT-integrated bins, we similarly calculate 
B^^dN/dy = Nj,^/{AtNEVT^y). 

We calculate Ae to correct for the geometric accep- 
tance of the detector and the inefficiencies of the MuTR 
and MuID, the track finding algorithm, and occupancy- 
related effects in the Au-|-Au environment. This is 
done by propagating PYTHiA-generated J/ -0 through the 
PHENIX GEANT-3 |16|] detector simulation, and embed- 
ding the resulting hits in real Au-|-Au events. The events 
are then reconstructed using the identical analysis as for 
real data, and the ratio is taken between reconstructed 
and embedded J/V'- The resulting Ae as a function of 
centrality is shown in Fig. [2] for both the North and South 
Muon Arms. The effect of the detector occupancy can 
be seen for more central events, as well as the higher oc- 
cupancy and resulting lower efficiency in the North Arm. 

The acceptance and efficiency as a function of px is 
relatively flat, with a modest 20% decrease from pT — 
to Pt ^^2.5 GeV/c^, then proceeds to rise again. The 
behavior is essentially the same across centralities, with 
the only difference being the absolute scale of Ae. 

We calculate the invariant yields separately from the 
two muon spectrometers, and then combine the values 
for the final results. We take the weighted average based 
on the statistical uncertainties and those systematic un- 
certainties which are uncorrelated between the two mea- 
surements. The final averaged result is assigned the un- 
correlated (reduced by the averaging) and correlated sys- 
tematic uncertainties. In addition, it was found that the 
invariant yields from the two spectrometers disagree be- 
yond their independent uncertainties, and a 5% system- 
atic was added to account for the difference. 



II. RESULTS 

The final J/ip invariant yields in 5%-wide centrality 
bins (integrated over all pr) are listed in Table HIl and as 
a function of pT in broader 20%-wide centrality bins in 
Table HTH 

The nuclear modification factor i?AA compares J /tp 
production in A -I- A with binary collision-scaled p + p 
reactions, and is calculated as: 



Raa 



1 dN^+^/dy 
(A^coii) dNP+P/dy 



(2) 



B 



d^N 



1 



N 



dp^dy 



2'KpTApT'^y Ae Nevt 



(1) 



where (A^cou) is the mean number of binary collisions 
in the centrality category of interest. The p + p results 
are from the combined analysis of data taken in 2006 
and 2008 as published in [ij]. The resulting Raa as a 
function of iVpart for J/^ from Au+Au collisions is shown 
as red circles in Fig. [3] 
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FIG. 1: (Color online) Top: muon spectrometer unlike-sign invariant mass distributions from same-event pairs (black, solid) and 
mixed-event pairs (red, dashed). Two example centralities and transverse momentum ranges are shown. Bottom: Corresponding 
unlike-sign invariant mass distributions after mixed-event combinatorial-background subtraction. Colored curves represent three 
different mass ranges for the signal extraction fit. 
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FIG. 2: (Color online) Acceptance x efficiency as a function 
of centrality for the North and South Arms. 



The systematic uncertainties are divided into three cat- 
egories: type A are the point-to-point uncorrelated sys- 
tematics, type B are correlated (or anti-correlated) point- 
to-point, and type C are 100% correlated (i.e. a common 
multiplicative factor) between all of the points. The error 
bars in Figs. [nHTU] represent the statistical and type A un- 
certainties added in quadrature, the boxes represent the 
type B uncertainties, while the type C systematics are 



included as text in the labels. The type A uncertain- 
ties are the RMS of the various mass fits as described 
above. The type B uncertainties on Raa are comprised 



TABLE II: J ftp invariant yields B^^dN/dy at forward rapid- 
ity (1.2 < \y\ < 2.2) vs. Au+Au collision centrality. The type 
C (global) uncertainty for all points is 10.7%. 



Centrality 


B^^dN/dy 


±stat 


±typo A 


-(-type B 


-type B 


scale 


(%) 






(Gov/ 








0-5 


1.25 


0.27 


0.11 


0.24 


0.12 


X10-* 


5-10 


9.12 


2.46 


0.13 


1.15 


0.81 


xlO"^ 


10-15 


9.37 


1.66 


0.71 


1.01 


0.95 


xlO"^ 


15-20 


9.16 


1.24 


0.48 


1.08 


1.08 


xlO^^ 


20-25 


7.11 


0.98 


0.29 


0.64 


0.64 


xlO"^ 


25-30 


7.85 


0.74 


0.29 


0.64 


0.64 


xlO^^ 


30-35 


6.14 


0.57 


0.28 


0.49 


0.49 


xlO^^ 


35-40 


5.43 


0.45 


0.14 


0.52 


0.52 


xlO^^ 


40-45 


5.07 


0.37 


0.15 


0.39 


0.39 


xlO"^ 


45-50 


3.49 


0.28 


0.07 


0.30 


0.30 


xlO"^ 


50-55 


2.76 


0.21 


0.07 


0.21 


0.21 


xlO^^ 


55-60 


2.85 


0.17 


0.04 


0.21 


0.21 


xlO~^ 


60-65 


1.64 


0.14 


0.03 


0.12 


0.12 


xlO"^ 


65-70 


1.17 


0.10 


0.01 


0.09 


0.09 


xlO^^ 


70-75 


9.49 


0.85 


0.11 


0.73 


0.73 


xlO"" 


75-80 


6.79 


0.69 


0.14 


0.51 


0.51 


xlO"*^ 


80-92 


3.43 


0.29 


0.04 


0.26 


0.26 


xlO"** 
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■ 2004Au+Au, |y|<0.35, global sys. =± 12% 
• 2007Au+Au, 1.2<|y|<2.2, global sys. =+ 9.2% 
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FIG. 3: (Color online) J/V' Raa as a function of A'^part- Error 
bars represent the statistical and uncorrelated systematic un- 
certainties, while the boxes represent the point-to-point corre- 
lated systematics. The global scale systematic uncertainties 
are quoted as text. The lower panel contains the ratio of 
forward rapidity to midrapidity for the points in the upper 
panel. 



TABLE III: J/iJj invariant yields 
ity (1.2 < \y\ 



""f'd^pi; forward rapid- 
< 2.2) vs. pT in four bins of Au+Au collision 
centrality. The type C (global) uncertainties are 10%, 10%, 
13%, and 19% for 0-20%, 20-40%, 40-60%, and 60-92% cen- 
trality, respectively. Bins in which the J ftp yield was less than 
the combined statistical and systematic uncertainties are cal- 
culated as 90% Confidence Level Upper Limits (CLUL). 





Pt 


R N 


stat 


type A 


+B 


-B 


SC3,1C 


{/o) 




dp^dy 




(Gcv/c- 






0—20 


0—1 


9.36 


1.41 


0.74 


1.01 


1.01 


V 1 n^^ 
X ±u 




1-2 


4.46 


0.66 


0.23 


0.40 


0.40 


xlO"** 




2-3 


1.37 


0.29 


0.08 


0.17 


0.17 


xlO""^ 




o — 41 


9 QQ 


1.12 


0.09 


0.27 


0.27 


X 10~ 




4— O 


2.05 


0.43 


0.14 


0.18 


0.18 


X 10 ^ 




5—6 




90% CLUL = 


3.27 xlO"** 






0— 1 




90% CLUL = 


2.00 xlO"** 




20—40 


0—1 


5.08 


0.54 


0.18 


0.67 


0.67 


X 10 




1—2 


2.78 


0.22 


0.09 


0.26 


0.26 


X 10^^ 




2-3 


1.11 


0.10 


0.02 


0.09 


0.09 


xlO""^ 




3-4 


2.76 


0.34 


0.11 


0.25 


0.25 


xlO"^ 




4-5 


7.47 


1.37 


0.31 


1.35 


1.35 


xlO"* 




5-6 


2.68 


0.61 


0.08 


0.31 


0.31 


xlO"** 




6-7 




90% CLUL = 


7.15 xlO"" 




40-60 


0-1 


3.19 


0.21 


0.06 


0.26 


0.26 


xlO"'' 




1-2 


1.49 


0.09 


0.03 


0.12 


0.12 


xlO"** 




2-3 


4.80 


0.38 


0.11 


0.39 


0.39 


xlO"^ 




3-4 


1.27 


0.13 


0.02 


0.11 


0.11 


xlO"" 




4-5 


3.86 


0.49 


0.02 


0.41 


0.41 


xlO"** 




5-6 


7.51 


1.69 


0.05 


2.04 


2.04 


xlO"" 




6-7 




90% CLUL = 


2.82 xlO"° 





of many sources, dominated by uncertainties in (Ncow), 
uncertainties in the matching of Monte Carlo and real 
detector performance, and differences in signal extrac- 
tion methods. The type C uncertainties are dominated 
by the normalization in the p + p invariant cross section 
measurement. Important systematics on the invariant 
yields are listed in Table HVl 

For comparison, in Fig. [3] we show our previously- 
published midrapidity J/^/i Raa values from data taken 
in 2004 The midrapidity measurement was made in 
the PHENIX central spectrometers via the J/ip dielec- 
tron decay. There is no PHENIX updated measurement 
at midrapidity from the 2007 data set due to significantly 
increased conversion backgrounds from this engineering 
run of the PHENIX Hadron Blind Detector The 
ratio of the new forward rapidity data to the previously- 
published midrapidity data, shown in the lower panel of 
Fig. El is in agreement with the previous result [Sj, where 
the latter led to speculation as to what mechanism could 
cause a narrower rapidity distribution in Au-|-Au than 
p + p collisions. 

We also calculate Raa as a function of pr, again using 
the published 2006 and 2008 p + p data 17]. Shown in 
Fig. 2] are the new results at forward rapidity along with 
the previously-published 2004 midrapidity results [5] . In 



0-1 


9.05 


0.57 


0.07 


0.73 


0.73 


xlO" 




1-2 


3.40 


0.22 


0.04 


0.27 


0.27 


xlO" 
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2-3 


9.19 


0.91 


0.16 


0.75 


0.75 


xlO" 
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3-4 


2.21 


0.35 


0.04 


0.20 


0.20 


xlO" 
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4-5 


8.13 


1.39 


0.01 


0.70 


0.70 


xlO" 
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2.31 


0.54 


0.00 


0.40 


0.40 


xlO" 
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7.69 xlO 
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TABLE IV: Systematic uncertainties on dN/dy for central 
and peripheral centrality categories. 



Source 


Central 


Peripheral 


Type 


Signal extraction 


9.8% 


1.3% 


A 


Acceptance 


3.4% 


2.2% 


B 


Input y, Pt distributions 


4% 


4% 


B 


Difference between mixed event/ 








like-sign background estimates 


1.5% 


0.6% 


B 


North/south arm agreement 


5% 


5% 


B 


MuID efficiency 


3.6% 


2.8% 


B 



some centrality bins for px > 5 GeV/c, we have no sig- 
nificant J/ip signal in Au-f Au and thus can only quote a 
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FIG. 4: (Color online) J/^p Raa as a function of px in four 
centrality bins. Error bars represent the statistical and un- 
correlated systematic uncertainties, while the boxes represent 
the point-to-point correlated systematic uncertainties. The 
global scale systematic uncertainties are quoted as text. 



90% confidence level upper limit on i?AA- 

As there has been much recent interest in whether 
the -Raa as a function of pr rises or falls, we have 
performed a simple linear fit to the -Raa at forward 
rapidity over the full pt range and obtain the follow- 
ing slope (m) values: m = -^0.011 ± 0.018 c/GeV (0- 
20% central), m = -1-0.065 ± 0.023 c/GeV (20-40% cen- 
tral), m = -f 0.034 ± 0.033 c/GeV (40-60% central), and 
TO = -0.037±0.053c/GeV (60-92% central). The quoted 
slope uncertainties are the quadrature sum of the statisti- 
cal and systematic uncertainties. Thus, only the 20-40% 
centrality indicates a statistically significant increase in 
-Raa with pT- 



III. MODEL COMPARISONS 

As previously mentioned, various theoretical models 
have been proposed to reconcile the J/ ^p suppression pat- 
tern previously published @ . Here we compare our mea- 
surements with three calculations for the centrality and 
rapidity dependence of the suppression. The first deals 
entirely with initial-state effects, while the other two in- 
corporate strong final-state effects. Then we compare our 
measurements with a simple cold nuclear matter effect 
calculation extrapolated to Au-|-Au collisions. Finally, 
we compare to a model calculation for the pt dependence 
of -Raa at forward rapidity. 

In addition to the models discussed here, there are 
many more models that only have a midrapidity predic- 



tion for J/^ production. Because this paper is focused on 
forward rapidity J/ip production, we have not included 
comparisons to those models. 



A. Gluon Saturation Model 

In the first model by Kharzeev et. al [l^, it is assumed 
that the nuclear wave functions in very high-energy nu- 
clear collisions can be described by the Color Glass Con- 
densate (CGC). The primary effect is the suppression 
of J/ip production and narrowing of the rapidity dis- 
tribution due to saturation of the gluon fields in heavy 
ion collisions relative to p + p collisions. In addition, 
the production mechanism is modified from p + p such 
that the multigluon exchange diagrams are enhanced. It 
should be noted that this model does not include any hot 
medium effects, but does have a free parameter for the 
overall normalization factor for the Au-|-Au production, 
which is fixed to match the midrapidity central collision 
J/V' suppression. Thus, the suppression trend with cen- 
trality and the relative suppression between mid and for- 
ward rapidity are predicted, but not the overall level of 
suppression. 

The resulting -Raa values calculated using this model 
are shown in Fig. [5l This model provides a reasonable 
description of the data and in particular matches the ob- 
served larger suppression at forward rapidity than mid 
rapdity in central events (-Raa'"'^/-Raa ~ 0.5). It is no- 
table that this ratio is essentially independent of central- 
ity in their calculation, whereas the experimental data 
shows the relative suppression approaching one in the 
most peripheral events. Additionally, the calculation at 
midrapidity actually indicates a significant enhancement 
(i.e. -Raa > 1) for peripheral events with A^pait < 50. 
This enhancement is related to a coherence effect of dou- 
ble gluon exchange. However, the coherence predicts an 
enhancement in c?-t-Au collisions and -R^Au at midrapid- 
ity, and no such enchancement is seen in the experimental 
data [l7| . 

Recently, the appropriate normalization factor for the 
above CGC calculation has been calculated [13, HH . Re- 
placing the normalization factor previously applied to 
match the theory to the magnitude of the observed 
midrapidity suppression, results in a predicted CGC sup- 
pression approximately a factor of two smaller than in 
the /auau/ data. This result suggests the importance of 
additional hot nuclear matter effects. 



B. Comover Interaction Model 

The second calculation comes from the Comover Inter- 
action Model (CIM) [l^,!!^. This calculation uses a rate 
equation that accounts for J/ip breakup due to interac- 
tions with a dense co-moving final-state medium. Addi- 
tionally, the contribution from interactions with the out- 
going nuclei is included. No assumption is made about 
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FIG. 5: (Color online) J/t/j Raa as a function of A'^part- Model 
calculation from [lH. Lower panel is the ratio of forward to 
midrapidity points and curves from the upper panel. 



the nature of the co-moving medium, i.e. whether it is 
partonic or hadronic, only that it can be represented by 
a comover density and comovei- J / cross section Cco, 
for which a value of 0.65 mb was found to match the 
NA50 data and then used for the projection to Au+Au 
at RHIC. The separate nuclear breakup cross section was 
taken to be tTbr = 4.5 mb. This value was taken from 
measurements in p+A collisions at y/s^ = 27.4/29.1 
GeV at the CERN-SPS. Under the assumption that (Tbr is 
energy-independent, the value from those measurements 
was used until recently as the cold nuclear matter refer- 
ence for heavy ion collisions at y^s^ = 17.2 GeV [2l|, as 
well as in [22| as the reference for ^s„„ = 200 GcV at 
RHIC. 

However, the effective breakup cross section has now 
been shown to decrease significantly with collision en- 
ergy [13]. In fact, a recent measurement in p+A colli- 
sions at ^s„„ =17.2 GeV Q yielded a value of tTbr — 
7.6 ± 0.7 ± 0.6 mb with no anti-shadowing correction, 
and CTbr — 9.3 ± 0.7 ± 0.7 mb with anti-shadowing cor- 
rected for. The value measured in d+Au collisions at 
RHIC [IHl is 2.8^i['4 mb (after shadowing is accounted 
for using EKS98 nuclear PDFs), somewhat smaller than 
the 4.5 mb used in The calculation, shown in Fig. [6] 
as the black, dot-dot-dashed curve, significantly over- 
estimated the suppression measured at midrapidity for 
Au-I-Au collisions at y/s^ = 200 GeV d, Q- The sup- 
pression is stronger than the SPS case mainly due to the 
larger comover density calculated for RHIC. 
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FIG. 6: (Color online) J/il^ Raa as a function of A^'part. 
Curves are calculations within the Co-mover Interaction 
Model (CIM). The black dot-dot-dashed curve is the older 
CIM calculation [2^ 1 pr edicted from SPS data. The remain- 
ing curves are from [23], where the dot-dashed curve is from 
shadowing alone, the dashed line also includes dissociation in 
the co-moving medium, and the solid line is the total effect 
after including J ftps from regeneration. The thin blue curves 
are calculated for midrapidity, while the thick red are for for- 
ward rapidity. The lower panel contains the ratio of forward 
rapidity to midrapidity for all points and curves in the upper 
panel where both are calculated. 



An updated calculation [23| was then released that re- 
placed the constant nuclear breakup cross section with a 
Bjorken-x-dependent function that goes to (Tbr — mb 
at ?; = 0, while the same CTco = 0.65 mb was used for 
the comover interactions as before. Additionally, a J/^p 
regeneration component was added that is normalized to 
the ratio of open charm production squared to J/ip pro- 
duction in p+p collisions. These new results are included 
in Fig. [51 The suppression from initial-state effects alone 
is much weaker at midrapidity than the previous calcula- 
tion, due to both the change in the nuclear absorption, as 
well as an updated parametrization of shadowing effects. 

The CNM effects (i.e. shadowing and nuclear absorp- 
tion) are much stronger at forward rapidity than midra- 
pidity, due in part to the assertion that nuclear absorp- 
tion is negligible at midrapidity. On the other hand, 
the effects of comover dissociation and regeneration are 
stronger at midrapidity. The combination of these three 
effects leads to predictions which are overall very similar 
at forward and midrapidity (as seen in the lower panel). 



10 




50 100 150 200 250 300 350 400 



FIG. 7: (Color online) J /ijj Raa as a function of A'^part- Model 
calculations by Zhao and Rapp from [SJ, [l^ are included for 
both rapidity bins, incorporating cold and hot nuclear matter 
suppression as well as coalescence of cc pairs. The various 
line styles represent the different contributions to the total 
as laid out in the legend, while the two thicknesses represent 
the two rapidity ranges (thin blue is midrapidity and thick 
red is forward rapidity). The lower panel contains the ratio 
of forward rapidity to midrapidity for all points and curves in 
the upper panel. 



C. QGP/Hadron Gas Model 

The third model we compare with is from Zhao and 
Rapp [H, [111, which incorporates both a quark-gluon 
plasma (QGP) phase and a hadronic gas (HG) phase. 
In this calculation, they include two different models for 
cold nuclear matter effects. In the first case, nuclear ab- 
sorption is calculated in the usual Glauber formalism, 
shadowing plus anti-shadowing are assumed to roughly 
cancel, such that the overall shadowing effects are encap- 
sulated in the breakup cross section tTbr, and px broad- 
ening is included via Gaussian smearing. In the second 
case, the cold nuclear matter effects are treated as in p3| . 
supplemented with the same px broadening model as the 
first case. 

The thermal dissociation is modeled via a Boltzmann 
transport equation for both QGP and HG phases. The 
QGP is assumed to be an isentropically expanding cylin- 
drical fireball. J/^/^-medium interactions are assumed to 
stop at a freeze-out temperature of 120 MeV. The final 
J/'f/' pt distribution is calculated by spatially integrat- 
ing the final phase-space distribution. The regeneration 



component assumes that the cc is thermally equilibrated 
with the medium when it coalesces into a J/^. Con- 
sequently, and the t/j pT distribution is governed by a 
blastwave equation for the transverse flow velocity. The 
normalization of this component is performed by plug- 
ging the initial charm densities into a rate equation with 
both gain and loss terms, and solving at the freeze-out 
time. 

The calculation results, using the second case for the 
cold nuclear matter effects, are shown in Fig. [3 along 
with the separate dissociation and regeneration compo- 
nents. Though the second cold nuclear matter case in- 
creases the suppression for more central events compared 
to the first scenario, the difference in the overall suppres- 
sion between the two scenarios is small. The qualita- 
tive trends of the calculation agree with the experimental 
data; however, the calculated suppression is very similar 
between forward and midrapidity, which is in disagree- 
ment with the data. 

It is noteworthy that the regeneration component is 
only slightly larger at midrapidity in this model than at 
forward rapidity. This is in contrast to other regenera- 
tion or recombination calculations that result in a signifi- 
cant narrowing of the J/ip rapiditydistribution in central 
Au-I- Au events (see for example [27l, f28| ) . In simple cal- 
culations, the J/ip recombination contribution scales as 
the square of the local charm density {{dNcc/dy)"^) and 
thus there is substantially less recombination at forward 
rapidity. This modeling also leads to predictions of signif- 
icantly larger recombination enhancements at the LHC 
where charm production is much larger. However, in this 
calculation 25] with a full space-momentum distribution 
of charm pairs, the probability of a charm quark from 
one cc pair recombining with an anti-charm quark from 
another cc pair is suppressed because they are typically 
spatially separated which is then maintained through col- 
lective flow. Thus, their recombination is dominated by 
the case where a cc is produced as a pair that would 
normally not form a J/ip, but due to scattering in the 
medium have a re-interaction and recombine. In this 
case, the regeneration contribution has a rapidity depen- 
dence similar to that of the directly produced J/ip- 



D. Shadowing/Nuclear Absorption/Initial-state 
Energy Loss Model 

In addition to the above three models, we use a frame- 
work for calculating just the cold nuclear matter effects 
and extrapolating them to Au-|-Au collisions. We be- 
gin with the prescription in [26| for d+Au collisions, 
which combinines effects of nuclear-modified parton dis- 
tributions functions (nPDFs) using the EPS09 parame- 
terization [29I with a rapidity-independent J / ip-micleoii 
breakup cross section tJbr, along with the possibility of 
initial-state parton energy loss. We have extended these 
calculations to the Au-|-Au case using the identical code. 
First, we include the variations of the EPS09 nPDF sets 
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FIG. 8: (Color online) J/-i/) Raa as a function of A'part. Shown in the left panel is a comparison with a calculation including 
cold nuclear matter effects (nPDF and atreakup) following the prescription in 17, '26^. The green (magenta) lines (light gray 
(dark gray) in B&W) are the result at midrapidity (forward rapidity) for all 31 EPS09 nPDF variations and the labeled value 
for atreakup ioT 0, 3, 6, and 9 mb. In the lower panel, one sees all 31 EPS09 variations times the various atreakup values resulting 
in less than a 10% difference between forward and midrapidity modification. The right panel shows the same calculation, but 
now including initial state parton energy loss with a quadratic length dependence. The chosen strength of the energy loss is 
that which most closely matches the d+ An data as detailed in paper [2^ . The lower panels in both cases contain the ratios of 
forward rapidity to midrapidity for all points and curves in the upper panels. 



and the breakup cross section crbr • Shown in the left panel 
of Fig. ISl is the projected cold nuclear matter effect from 
these two contributions. The top green band shows the 
results for J/^ at midrapidity from all 31 EPS09 nPDF 
variations with a crbr = mb. The lower green bands are 
in steps of (Tbr — 3,6,9 mb. The magenta bands show 
results for J/ip at forward rapidity. The centrality de- 
pendence is not reproduced for any value of (Tbr at either 
rapidity, though a value of CTbr greater than 6 mb (9 mb) is 
required to approach the data at mid (forward) rapidity. 
In the lower panel we show all 31 EPS09 nPDFs times 4 
(Tbr values spanning the range a^r = 0-9 mb for the ratio 
of the forward to midrapidity suppression. No combina- 
tion of these two effects reproduces the modification in 
the rapidity shape for mid to central Au-|-Au collisions. 
One reason for the modest rapidity dependence is that 
at forward rapidity the Jf-tjj production results from one 
low- a; gluon (in the nPDF shadowing regime) and one 
high-x gluon (in the nPDF ant i- shadowing regime) and 
the two effects largely cancel. 

As discussed in [26|, one can attempt to improve the 
cold nuclear matter calculation agreement with the d-\-An 
data by including a parameterization of initial-state par- 
ton energy loss. In the right panel of Fig. [51 we again plot 
all 31 EPS09 nPDF parameterizations, 4 (Tbr values and 



the quadratic length dependent initial-state parton en- 
ergy loss that best matched the d+An data. The initial- 
state parton energy loss has a minimal effect over this 
rapidity range, which is not unexpected since the effect 
only becomes significant in d+Au for rapidity y > 1.8. 

We note that the cold nuclear matter calculation does 
not give a full description of the d+Au data, adding some 
uncertainty to its use in Au-|-Au collisions. Nevertheless, 
it is informative that the calculation clearly fails to si- 
multaneously explain the Au-|-Au data at forward and 
midrapidity. We observe J/ip suppression beyond that 
expected from the cold nuclear matter effects included in 
this calculation with the choice of a reasonable value of 
3.0-3.5 mb for (Tbr at RHIC. 



E. pt Dependence of the Suppression 

Most of the above calculations do not include predic- 
tions for the J/ip suppression as a function of trans- 
verse momentum. However, the calculation of Zhao and 
Rapp 25] provides nuclear modification factors at both 
mid- and forward rapidity as a function of px- Shown 
in Fig. [S] are the results compared with our experimental 
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FIG. 9: (Color online) J/^p Raa as a function of pr in 
four centrality bins. Model calculations by Zhao and Rapp 
from [2^, [2^ are included for the forward rapidity bin, incor- 
porating cold and hot nuclear matter suppression as well as 
coalescence of cc pairs. 



data. Their calculations indicate a moderate rise in Raa 
versus px predominantly due to the Cronin effect [s^ . In 
fact, in other recombination models the enhancement is 
limited to low px p7l | , and in this calculation the recom- 
bination contribution drops off beyond pt ^ 3 GeV/c. 
At low pt in the most central bin, the suppression from 
this calculation is too weak by up to a factor of two. 

More recently, Zhao and Rapp [13, HH have modified 
this calculation to include feed-down from B-mesons and 
a reduced suppression at higher px due to the longer 
formation time of the preresonance state to the J/ from 
time dilation. These contributions serve to increase Raa 
at higher px compared to the previous calculation, and 
are compared to the forward rapidity data in Fig. 1101 
The current lack of statistics in the data at px > 5 GeV/ c 
precludes a confirmation of this effect. 

It should also be noted that in the new calculation by 
Zhao and Rapp, cold nuclear matter effects are handled 
differently than in the earlier calculation. An effective 
absorption cross section of 3.5 (5.5) mb at y — (1.7) 
is used to account for the combined effects of shadowing 
and breakup. These effective cross sections are obtained 
from comparison with recent PHENIX d+Au data. They 
argue that the larger effective breakup cross section at 
forward rapidity is most likely associated with shadow- 
ing effects, and so refiects a suppression of the number 
of charm pairs relative to midrapidity. Therefore, the 
additional effective absorption at forward rapidity is as- 
sociated with a reduction in the open charm yield as well 
as the J/ip yield, thus also reducing the J/ip regeneration 
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FIG. 10: (Color online) J/ip Raa as a function of pr in the 
0-20% centrality bin at forward rapidity. Overlaid are a more 
recent calculation of Zhao and Rapp [sdi, (sH) . 



contribution. As a result of the ad hoc use of a larger ef- 
fective absorption cross section at forward rapidity, this 
new calculation produces a forward to midrapidity Raa 
ratio of about 0.7 in central collisions, which is in bet- 
ter agreement with the data. Again, this is entirely due 
to the cold nuclear matter effects, and in fact the hot 
nuclear matter suppression in this calculation is almost 
identical between forward and midrapidity. 

A second model of interest to the px dependence of 
J/ip production is the so-called Hot Wind model [33| . 
This model predicts a decrease in J/ip Raa at higher 
Px in semi-central events, based on a modification of the 
screening length due to the relative velocity between the 
J/^ and the medium. However, there is no quantitative 
calculation available at forward rapidities, and there is 
no evidence for such an effect in the px range covered by 
the present data. 



IV. SUMMARY AND CONCLUSIONS 

We have presented new and more precise measure- 
ments of J/ip nuclear modification at forward rapidity 
in Au-I-Au coUisions at ^s„„ — 200 GeV. The results 
confirm our earlier published findings of a larger sup- 
pression at forward compared with midrapidity. This, 
combined with the similar suppression of J/ ip at midra- 
pidity between RHIC and lower energy measurements, 
remains an outstanding puzzle in terms of a full theoret- 
ical description. 

Due to the lack of a comprehensive and consistent un- 
derstanding of the numerous cold nuclear matter effects 
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and their extrapolation to nucleus-nucleus collisions, ex- 
tracting a quantitative measurement of the hot nuclear 
matter effects is not possible at the present time. How- 
ever, it is clear that the observed suppression is at- 
tributable to hot nuclear matter effects, as the suppres- 
sion in Au4-Au collisions is larger than predicted by the 



current models of CNM effects. It will be useful for all 
calculations to include the transverse momentum depen- 
dence for future comparisons, especially as the experi- 
mental uncertainties, particularly at high pT, will only 
improve in the future. 
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